{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Ray tracing: Cython with structs"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "%load_ext Cython"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%cython\n",
    "cimport cython\n",
    "import numpy as np\n",
    "cimport numpy as np\n",
    "DBL = np.double\n",
    "ctypedef np.double_t DBL_C\n",
    "from libc.math cimport sqrt\n",
    "\n",
    "cdef int w, h\n",
    "\n",
    "cdef struct Vec3:\n",
    "    double x, y, z\n",
    "\n",
    "cdef Vec3 vec3(double x, double y, double z):\n",
    "    cdef Vec3 v\n",
    "    v.x = x\n",
    "    v.y = y\n",
    "    v.z = z\n",
    "    return v\n",
    "\n",
    "cdef double dot(Vec3 x, Vec3 y):\n",
    "    return x.x * y.x + x.y * y.y + x.z * y.z\n",
    "\n",
    "cdef Vec3 normalize(Vec3 x):\n",
    "    cdef double n\n",
    "    n = sqrt(x.x * x.x + x.y * x.y + x.z * x.z)\n",
    "    return vec3(x.x / n, x.y / n, x.z / n)\n",
    "\n",
    "cdef double max(double x, double y):\n",
    "    return x if x > y else y\n",
    "\n",
    "cdef double min(double x, double y):\n",
    "    return x if x < y else y\n",
    "\n",
    "cdef double clip_(double x, double m, double M):\n",
    "    return min(max(x, m), M)\n",
    "\n",
    "cdef Vec3 clip(Vec3 x, double m, double M):\n",
    "    return vec3(clip_(x.x, m, M), clip_(x.y, m, M), clip_(x.z, m, M),)\n",
    "\n",
    "cdef Vec3 add(Vec3 x, Vec3 y):\n",
    "    return vec3(x.x + y.x, x.y + y.y, x.z + y.z)\n",
    "\n",
    "cdef Vec3 subtract(Vec3 x, Vec3 y):\n",
    "    return vec3(x.x - y.x, x.y - y.y, x.z - y.z)\n",
    "\n",
    "cdef Vec3 minus(Vec3 x):\n",
    "    return vec3(-x.x, -x.y, -x.z)\n",
    "\n",
    "cdef Vec3 multiply(Vec3 x, Vec3 y):\n",
    "    return vec3(x.x * y.x, x.y * y.y, x.z * y.z)\n",
    "\n",
    "cdef Vec3 multiply_s(Vec3 x, double c):\n",
    "    return vec3(x.x * c, x.y * c, x.z * c)\n",
    "\n",
    "cdef double intersect_sphere(Vec3 O,\n",
    "                      Vec3 D,\n",
    "                      Vec3 S,\n",
    "                      double R):\n",
    "    # Return the distance from O to the intersection of the ray (O, D) with the\n",
    "    # sphere (S, R), or +inf if there is no intersection.\n",
    "    # O and S are 3D points, D (direction) is a normalized vector, R is a scalar.\n",
    "    cdef double a, b, c, disc, distSqrt, q, t0, t1\n",
    "    cdef Vec3 OS\n",
    "\n",
    "    a = dot(D, D)\n",
    "    OS = subtract(O, S)\n",
    "    b = 2 * dot(D, OS)\n",
    "    c = dot(OS, OS) - R * R\n",
    "    disc = b * b - 4 * a * c\n",
    "    if disc > 0:\n",
    "        distSqrt = sqrt(disc)\n",
    "        q = (-b - distSqrt) / 2.0 if b < 0 else (-b + distSqrt) / 2.0\n",
    "        t0 = q / a\n",
    "        t1 = c / q\n",
    "        t0, t1 = min(t0, t1), max(t0, t1)\n",
    "        if t1 >= 0:\n",
    "            return t1 if t0 < 0 else t0\n",
    "    return 1000000\n",
    "\n",
    "cdef Vec3 trace_ray(Vec3 O, Vec3 D,):\n",
    "\n",
    "    cdef double t, radius, diffuse, specular_k, specular_c, DF, SP\n",
    "    cdef Vec3 M, N, L, toL, toO, col_ray, \\\n",
    "        position, color, color_light, ambient\n",
    "\n",
    "    # Sphere properties.\n",
    "    position = vec3(0., 0., 1.)\n",
    "    radius = 1.\n",
    "    color = vec3(0., 0., 1.)\n",
    "    diffuse = 1.\n",
    "    specular_c = 1.\n",
    "    specular_k = 50.\n",
    "\n",
    "    # Light position and color.\n",
    "    L = vec3(5., 5., -10.)\n",
    "    color_light = vec3(1., 1., 1.)\n",
    "    ambient = vec3(.05, .05, .05)\n",
    "\n",
    "    # Find first point of intersection with the scene.\n",
    "    t = intersect_sphere(O, D, position, radius)\n",
    "    # Return None if the ray does not intersect any object.\n",
    "    if t == 1000000:\n",
    "        col_ray.x = 1000000\n",
    "        return col_ray\n",
    "    # Find the point of intersection on the object.\n",
    "    M = vec3(O.x + D.x * t, O.y + D.y * t, O.z + D.z * t)\n",
    "    N = normalize(subtract(M, position))\n",
    "    toL = normalize(subtract(L, M))\n",
    "    toO = normalize(subtract(O, M))\n",
    "    DF = diffuse * max(dot(N, toL), 0)\n",
    "    SP = specular_c * max(dot(N, normalize(add(toL, toO))), 0) ** specular_k\n",
    "\n",
    "    return add(ambient, add(multiply_s(color, DF), multiply_s(color_light, SP)))\n",
    "\n",
    "def run(int w, int h):\n",
    "    cdef DBL_C[:,:,:] img = np.zeros((h, w, 3))\n",
    "    cdef Vec3 img_\n",
    "    cdef int i, j\n",
    "    cdef double x, y\n",
    "    cdef Vec3 O, Q, D, col_ray\n",
    "    cdef double w_ = float(w)\n",
    "    cdef double h_ = float(h)\n",
    "\n",
    "    col_ray = vec3(0., 0., 0.)\n",
    "\n",
    "    # Camera.\n",
    "    O = vec3(0., 0., -1.)  # Position.\n",
    "\n",
    "    # Loop through all pixels.\n",
    "    for i in range(w):\n",
    "        Q = vec3(0., 0., 0.)\n",
    "        for j in range(h):\n",
    "            x = -1. + 2*(i)/w_\n",
    "            y = -1. + 2*(j)/h_\n",
    "            Q.x = x\n",
    "            Q.y = y\n",
    "            col_ray = trace_ray(O, normalize(subtract(Q, O)))\n",
    "            if col_ray.x == 1000000:\n",
    "                continue\n",
    "            img_ = clip(col_ray, 0., 1.)\n",
    "            img[h - j - 1, i, 0] = img_.x\n",
    "            img[h - j - 1, i, 1] = img_.y\n",
    "            img[h - j - 1, i, 2] = img_.z\n",
    "    return img"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "w, h = 400, 400"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "podoc": {
     "output_text": "<matplotlib.figure.Figure at 0xcbb4a90>"
    }
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAc8AAAHPCAYAAAA1eFErAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAAWJQAAFiUBSVIk8AAAIABJREFUeJztvUusLF2W3/XPc+9X1WVHd5VM+90Dy0IwoEctwRCQQDzE\nACEMksXDAxjY2AJjJARYQhaSASFhDLKxBzAwD1kCI8QA8RBIwBAkj8wAywJLbuM22E1VV3RVfd+9\n9xwGcfeNFSvW2o/IzHMyz/n9pFRE7NjxynMyf7nW3rHj9PT0JAAAAOjn4aVPAAAA4N5AngAAAIMg\nTwAAgEGQJwAAwCDIEwAAYBDkCQAAMAjyBAAAGAR5AgAADII8AQAABkGeAAAAg7wfqfyDH/yAsfwA\nAOBV8t3vfvfUW5fIEwAAYBDkCQAAMAjyBAAAGAR5AgAADII8AQAABhnqbZvxve997xK7AQAAuDrf\n//73z94HkScAAMAgyBMAAGAQ5AkAADAI8gQAABgEeQIAAAyCPAEAAAZBngAAAIMgTwAAgEGQJwAA\nwCDIEwAAYBDkCQAAMAjyBAAAGAR5AgAADII8AQAABkGeAAAAgyBPAACAQZAnAADAIMgTAABgEOQJ\nAAAwCPIEAAAYBHkCAAAMgjwBAAAGQZ4AAACDIE8AAIBBkCcAAMAgyBMAAGAQ5AkAADAI8gQAABgE\neQIAAAyCPAEAAAZBngAAAIMgTwAAgEGQJwAAwCDIEwAAYBDkCQAAMAjyBAAAGAR5AgAADII8AQAA\nBkGeAAAAgyBPAACAQZAnAADAIMgTAABgEOQJAAAwCPIEAAAYBHkCAAAMgjwBAAAGQZ4AAACDIE8A\nAIBBkCcAAMAgyBMAAGAQ5AkAADAI8gQAABgEeQIAAAyCPAEAAAZBngAAAIMgTwAAgEGQJwAAwCDI\nEwAAYBDkCQAAMAjyBAAAGAR5AgAADII8AQAABkGeAAAAgyBPAACAQZAnAADAIMgTAABgEOQJAAAw\nCPIEAAAYBHkCAAAMgjwBAAAGQZ4AAACDIE8AAIBBkCcAAMAgyBMAAGAQ5AkAADAI8gQAABgEeQIA\nAAyCPAEAAAZBngAAAIMgTwAAgEGQJwAAwCDIEwAAYBDkCQAAMAjyBAAAGAR5AgAADII8AQAABkGe\nAAAAg7x/6RMAeGtM0/+o5aP3YF4nM7UvBVM7//TlNc+/cN0TB4AvIE+ACzJN/42kryS90yLDMi2y\ns7Isrx5pZiJdhTpNfz44oyfz+qR5/vnD1wYAK6RtAQAABiHyBDjANP3Xn+fea03BlijTpmHtvLSP\nMG3UGa1XsM6XKZjPzvv/+jz3JOnx8/yj5vlvam4LACvIE6DBNP0ZSd/WXpTStt3St11aMWbyzF5S\nLs2aQHuWn3bL0/SXtAr1UdIHzfPfLACIQZ4AAdP0X2qV5Ttto8rykvqkWYssa+I81vY5Go2uPJnp\ntzVNv6i1rfS3DewH4PWDPOHNM03/6ee5n9LS2acI06difYQptYXZijyjsqORZy3iHJGotBXpk6bp\nr6iIdOGD5vm3D+4T4PWAPOFNMk3/kaTvaI0upTXCzITpRSnVpemXpbooR9s6j0Sf50jUTr+lafol\nLTL9QGQKbw7kCW+KafrT2kaX5SXVo0x/i4mNJLNIsyVJv3zptG1vm2hWvm8b3c4/aXkfv6Vp+qta\n0ru/Jdk3wOsCecKrZ5r+lKRvKZZm1H7Z6gBkI89MlL0RZkuuStarMm2ti5azsoIXaSmzEemjpumv\naYlGP2qef2tlfwD3DfKEV8tWmu8US7PMS9sIs6cTkIKyUTleWqBR2WgKd0SiVp7vzPQrI1Jpnn9j\nZZ8A9wfyhFfDNP37kn6tFllKdWkWcdoUbE+0GaVte19ZfV+uYL5WZqe1dS1xtmT6VKljJfqgbVpX\nmqZf1iLSrzXPP5fsA+B+QJ5w90zTn9DS+aekZX0HoFqvWZu2bUmz1gGoFXnK1e+NQP20Z52S+Zo8\nL5nW9encMn0n6f3niPSD5vk3V/YHcNsgT7hbpumPay/NqANQT0cgmWUvxJFbTy4lUCXTS8jziDhr\n5T4i9eIs8yUiLRL9ZUnfaJ5/U7JfgNsFecLdMU3/rqRfo6U900vTyrMWadZuO+m5Z7OWhh2Vpy9T\nZZtoXkHd2rQ237Ps8eu9NO28leijpHdIFO6Sh3YVAAAAsBB5wl0xTX9S+W0n/tYTn649BfM+Bdsb\ndY5Gnkrq+7JoWcH6rMyvj6Z+fmRdD6V+LZ0bp3Hn+dcNHgvgZUCecBdM0x/TMjh71L4ZtW1K/eL0\nadsecfbKU6aerx9tr2TbnpTtiECvKc9oOytS36HoQUsK90HT9CuSPkoSIoWbBnnCzTJNf1TLeLPS\nVprR2LM+ypT6xek7/1xi0AMF61XZ1tdXMj03+szKsjrZcoSVY23dSfsotLzvjyp/j0WkHzTPf0PH\nsQGeF+QJN8c0/Zta7tcsjwGT+qLNSJ4tcY5Gka1lBeWt9bVlX19JnZGo88h8tJxRq1fE2YpC1+Vp\n+oGW+0N/Q+fxAa4P8oSbYpr+iBZp2p60Uty22XriiZfkJXrO1ujZphVtSvE2fr4WoUZllxLoJfD7\nt9GnF2v5G73TNH1f8/y9K5wPwDjIE26CafrXtUSbUWcgKe4UFEmzfNlKeZQZSVTJ/EgkmkmzFiG2\njt+zj1ZZTaCtfTwH9tp8JFr+XqVN9IdaotCffcbzA9jDrSrw4kzTvyVp0rZDkI80I3FGUWRve2VL\ngKqs641KR6PPaD5bNxqN1uTdOv5zYv+O9lX+/u8lfedzKhfg5SDyhBdjmv6Q9tI8Z8CD3qjQr1Ol\n/gi1bVtRn93ez0f1srJo/Wj9W6G8B3GnoiUK/Ynm+de/3CnCmwV5woswTX9Y69B6NXFGt51kqVo7\nL+1l1hM1nhu1KtlXVF5b9vV9+VGB9mx7a0QSXf/W0/R9SaI9FJ4V5AnPyjT9G5/nfKeg7Dmbvv3S\nv7JyJeW16PMSAhmJPHuiy9b2Ud3W+bWOdatYifr/AWmafqh5/ukXOzt4WyBPeBam6Q9K+mktwpT2\n49LWxp5VUH4kivRltl5UfzTy9NTKs+VLRqCt6b2S/y2madbSoYh7Q+G60GEIAABgEOQJV2ea/lWt\nUadt42z1ph1N1fZEiHLLl4zCeiPNqKwn4myVRetr0ealr/+5Kf8L9v/nKy29cb//kicGbwDkCVdl\n6VFbHh8WdQyK0rSj6dNeafakdo8exxKV19LF0fKRlG2PDF9L6rZQ/mf87Szf/jy8H8B1QJ5wNRZx\n/pTqT0GpPfmkJwrtaYe8dtQ5IlEpFmC2nB2vZ/1ra+tsYSPR95K+9bkNFODyIE+4Cqs4W48OywZv\nH40upVySCrbxdc6JNKPj+fJaWbbPVqRqy0bE+Jol6iPRrzRNv/qypwSvEnrbwkWZpj8g6WcUjxbU\n88STnmiyp56SZb+vI/RIdLRsNBrtPfe3EnV67P+TPgv0J/TChYtB5AkXY5r+BW07BrXGoz3pWDpW\nybQWZV46Vev3G5Vn9c89XlYWvT+m1kk6nU6fp/XX68CmcUtHov/vZU8JXg3IEy7CIs5J1xHn0bRt\nVmbLL5GurR27xeWiz0sK8fXItfyfPWjJgvwUPXHhIpC2hbNYpCktT0Qpt6G02jalrThbY9TKLbci\nrSz6Ojf67JVobTkrGz2PmsRq0fll8Md+ip5/fTOctPwvlpP+9peB5ef5uy91UnDnIE84zDT981rG\np5X6xenHnd22TfVFl35aizxr5UcZ3d+IcPPoM5bl5cV4BHtutytS+7+3ME0/QKBwCNK2cIhp+v1a\n798sQ+31iNNKsZW2VTAfTf28LcukcolUbe24R7ZLan5JlbauN9r/88v1ttO7NoVbbmfh8WYwDvKE\nYabpn9X6RBQ/sPuR9s1CFm3W1rX2oc51GUdStVFZaz/rOt+pZ0yUt2Ws2xSp/d9EoHAM0rYwxDT9\nPq0RZxGmtL23rva4sJYYo/RtVD/bh+eIMCOu8e3fk5atbVOLqrfTc+R1qTTsbaV2bdywfA1O0/d5\nrBl0Q+QJAAAwCPKEbqbp92iJOrOHV4+2afa2c9baP/28OsqP0pd2HT1+nNaM0r7Vvez2Z+/rPIdr\n3LJyO6lcmzH5NveBQjfIE7qYpn9a6+0otfbNc+7lVDCvgXm7/bW+mS+b/t0KpHffWfvoSadT9N51\n7PGgzM4V6e0ItPwA/Lam6a+/8PnAPUCbJ3Ty04rHqT36sOpIcq12zaxOtuwZEXdrPyPlQc0vVUei\n2ZPZtv0+nA6YqXeTqM3yaJtm2e5l20EfJD2pDKQA0AJ5QpOld23pIBRFndL1U7N+2c+3JDQikuuF\nQ2055RVq0tyLMt7PkUivJcqozpFBFF5eouV/+b2m6Vc1z7/2pU4E7gDkCVWm6fdqGeS9iDMb9CDr\nXSszH0msJktfrzUfbXNUhJcSqI0Wj+07E2Mu0752zlodK7CsXlanJtvblqj9f5amadY8T899EnAn\nIE9Imabfrb04o+dvSntpZrKUqV9L0x5J20b7PYfL7Gc82ozkZ6NMv80pTdEekejTU3u7qE6RXU2k\nty/R8n+8pHCn6Vc0zz/znCcAdwLyhJBp+qe0faxYdA+nl+do2laVMjtVsJyJ87mo/Sj4PDd0Oj5C\n3W+8CNKnateySKC1c8jWHZGWl2lLpKMSfRmBSsvzQLn/E/YgT9gxTf+41hGEijDtrSi1sWprMvX0\nRJ61utn6l2c02qyNW+ulaKUZydZWj/bbK/Wa4GpCsyJtCbNXjM8v0PJ//U7LCERLD1yeBwoF5AkB\nUc/a1pNPjqRqFSxH27S4pjjH933k8V/rsXwEuq7LpLmXa3u+ds69qdZRofloNFuunePLRKDLLSwA\nFuQJX1ju5ZRycY6kYVuR50i61m9XW35Z6k8+8XWzNKtPx+73H4nU765VnhGtr6Vlyzovtqzt1G9r\nj9sjULuP61LkWW5hkabph5rnn36Og8ONgzxBkjRNv0vrr+vodpRTsBxFngrWKahjyxSUR2W3JUrP\nSMQ5lm5do04/tXVjybb23UckxjL1Aq0J1e8ji0KjY1qeLwq1ApWWTkS0gQIjDAEAAAxD5Amapn9U\ny6gq5d/BR54ntW9DaUWbvVFkT+q2tY9L8pQcp0QiWRtlzj6Nun1Pt9GlgrJ1XRZZjkaetfSqn7fL\ndjoSbfqoM9quta/njz6l5e/+labp/9U8//rnODjcKMgTtO0gJI2PT1voaeeMlu+7rbNPnLFkezr3\n2HbOVZ7rukyUrXW9RBLN5GlFmPW69fuOZFq4DYF6eb4TQ/gB8nzjTNM/oW0HIandu9a3b9qpJxJs\ni9sRY4vR2z5abZxxm+ZaZiNO+9ru47IC9aIs85E4ozL/PmQdi2q9cVsC9ed3ecqbZQdQ+IHm+bvX\nPCjcMMjzDTNN/4jWgRCiQd5b0nyOqDMqfy65ZinbscitfdvJPhW79rrdR5yRNLPlHon6c22laMt8\nTZi9Mo2IpNnbE/f5otAi0L+mef7Zax8UbhDk+aYpz+a0kaa0ym9UoJ4jUeftczrlUt3XzcvidspI\npiXqrEvy4cHva1+ntw200JOmbb1q9aPjZW2hz9fGWcN+Nh7E/Z9vF+T5Rpmm36nlSSnRE1JqkaeC\neQV1LCPSvBXZ9gsyo1+cPrJsy9MKsyZVu96vq51noVeej4/5sq3/+JgfJ0vXjgj0+aLPItB33Lry\nRkGeb5Bp+ge1DvjuB0KQ2tGmXJmntc7P309EOt7GuV/et2+uy/61ltfFWebL/r1Uswh0ez5rWa19\n05YVST48bJeLxIosa1Ho42MsPd+Z6HYEauVZ0rf/j+b5N1zzoHBjIM83yaQ1XetvSZHq6dqaGD23\nEkWOEkedi1zaEWldnD6q1JflXJ65LP201PflNYG2fhC00rBFmI+Py/GsDMu+Sx27bRFrEW927NsX\nKOnbtwjyfGNM0z+sfTtnFF3Wok5fJ1pWsE2LWxBsTY5HxbkdfzaKPNfo0spzqZeJ0k6tPO1yS6Lb\n83FX2+hNa6VZpkWcVqJl30WAXqRZKteex+0J1P6/l/TtX2fg+DcE8nxzlMEQStRpfz1n8lSwrka2\n/tYj0VyOPenatY6NLrP9nb5Mo9fDQ12eXpDlVepHYo0EWs47u74sWoyk6V/+GGVbO1+i1ZZEjwr0\numzlufwohbcC8nxDTNPv0D7qbEWeCqZy9SJ6RXkLEi3fvjVxjnUgitoUy/y+LXOdLtI8fZFd2cYK\nsvaS4nIv0iPyjCLG7PXp0zovrWVWnPYYNuWbUROoPVf/d3ie6HOR6DT9sub5113rgHBDMLYtAADA\nIESeb4pya4p9sPUl2zajslunJ6Icb+v05dt2zriN00edtcjz3bt4WcrrROnbsv8o+sw6B0nblK2N\nMsvyw8NaLq1RZ5nal0/hHok+W3+b60efJYvDV+pbgb/0G2HpKBQ9asynbguRJHtk+BzCPP8ezN79\n9AyIEHcS8tOTWe/F6SXaTsMWMfqptC+z81EHo6jzkJTfliLFKdoyLeL07ZxeoBFl+zIfkYnwZdo/\nt2lbUrdvB+T5Bpimv09L1PleuTCzeZmy10K9jXNbb+y6Y5H6ezn34lyizodmG6YVY3nZZWlfFkWh\nvk01a/f0UaeNEiN5lqizvD592h7DLmfYW1isQO2tL3bZCjMT6PXFagX6XtP0S5Kkef5N1zwovCDI\n803wa7QVZ5SujdK2WfQZbSNXp6fsHI5En73SXOqeat/wZU+NdG0R53ZfdtCDXJxWnpE0o5e0L4si\n1EyeXk7R/ZzSPtos4vz0Sfr4cX8N/hhForX30N7qUhNlrcPQ83DS+v9YPmPfeqmTgWcCeb5ipunv\n/TyXDYgwIsCsXla/Vd5DS5AjEeTIufSI039Tbx8Ttkz3+1jv34xvTfHtm1Iuzvfv98ulvl/vI9Ga\n2L5cYSBP33s2epWI8+PHfWq41avXzpf3ImsLzdo9nz/69J+lcuuKNE2/RPT5SkGer5rvfJ6WqPOk\n8cjTTm+V1jdivzTb9Zc6Xoy1dk874IFv31xeD2mnIGkrPivEMm+n0r48i0TLOUW3r0jb9k0vTx9x\nfvy4jTqLOKNjRO/X5h02f04/0EKtrpXjbQhUIgJ9vSDPV8o0/d1ab9qu9bCVm5cpU1Be41zJRpHm\naNR49LitY7TqLOduv+C3X/anz2VRxHnaidNHnpks7avINiov9X0Kt5yn70AkxVGnlad9vXu3T9fa\na/jwofLW+ncykJtP4W7e2Upq92Ww8nzQNP0VzfNvfskTgiuAPF8t39H6521FnYWeiLNXYiPR3ktF\ntrEQ9xHO9hyzKHO7bGWpnVRtBGqjMht1Rm2YkTi/+mqdl/by9BK1EaiUR57SXp6lrfLxcRttRvuO\nevRG+7fHsUMAlnK/bbYPz8tHn9LymWPkodcI8nyFTNPfruUDu/76badrs8jTLtci1UuRyfSSEeh1\n2kpjQa6SXMqj4fjiiM1HnlagVpp2Ku2FmgnU9+bN5FbStv6+zSLO8vrwYd8hqSdNW45Tpnb+4WEr\nULsvP8JQKXuZ6POk9X/Xy/OdpukXNc8/99wnBVcEeb5KJq33dErtiLOQfcNdSlgvvZ+2DJcv4rxe\nX9T5ZWlXvm3zlB4eThvZRG2eXppemOVVIk9fVhOoNCbPEnl6eRZxfviwl39EFnG+e7ddZ6dS/MSW\nl326yuZo2v/vlIwPT115bSDPV4kdv1aKI08FyzLl0byn1UZ4juRa++49B8tlI83NnjdR5/ZWFDtA\n/FaacQ/bKG3r2y69OMtLissjgfrevD496ts8vTxt1OnFGaWqN+90Is8iULtuG8XHkvf7bkmzR7zH\n8Z+xtfctvB6Q5ytjmv4ebXvXSsfStFnZKLU07KlSp1e+tW+/0e3r9Vupx2y9Tzm2Ura2zVDa97b1\nEehXX0nf+tZWnmU5E2g0qMKIPEt754cPqzy9OHuEWaLaqM2zpGwzidp9tv4uLxN9luVFoNP0VzXP\nv/G5TgSuTCWxAgAAABFEnq+O0lGovKQ44owi0N507TU4J/o8ejxt9p9HlvmACH67df7kIk6btrXD\n8q2RVC1tW4s6bfQprfM90ae0ve9zJPKM0rVZ2+nm3QyiztKb10fcPiK155gN12fnWxHndXve2vmT\nSN2+LpDnK2Ka/g5tB3/PRNmTsr0VLn2f52hP26huXeqrPJb33A6S4Ns+rTQjgUr1nrZeotIqzpZA\nfdrWi89Kq/S0/fhxWWflaXvv1lK2/n7R6GHaZVla5WkfrO1fPePaPj82dbtt+5ymv6x5/q0vc1pw\nUZDnq+Lb2g+GIKkronyOqNNLp7Xst1Vlfeu4hV5pbuvtI6ltvSjSikTie9vWOgzVBknw0iyilJbp\nt7+9laeVqI0+pe3xInnaezxt5Fl62LbaTP2+vDDfvdsvS9totNX2WY5h3+fbEav98cpX7muBv+Sr\noqRsWynaKI17hF75jaRfR3ra9lI/9vLYsb66S/183n+px2nbNWWbRZ2tyDNL3Ur7yDOKQK08o/Fu\nC16eNvKMBkSIola/n+gZoP4lrfuNok/7nmf3e/r1tb/nZcSaRZzlRa/b1wTyfCVM09+peOB3BfPX\npggwa8e05/HcbZ3+POy5HN1PvL0V57637VYGI4MkRPd82jZPG23aKNSL1g7n50UuxZGnl2ckXd9W\nKu2jVzs2brkm+zgzaZ23A8NHadvr3nZyKVaJTtP/rXn+LS99QnAmyPPVYDsKZVFlK2X7HGKN6JHs\npY5T6E3fevbbZSnbZbp9b6OUbSbRaGD42hB9/laVnuhTyofUizoLRe2kUbRqpSut+/jqq/2YuPZV\nBFqu2wq01e7ZStVeP3Xbs3Oiz9cC8nwFTNPfqu29nVmaVsn65+Boe+e5EeJRYa51t2nC+DyzVG5Z\nrqV3a+lbqS/ytGnYqC00EqiXp+1xW/Cp1hJ5lh62/lqj1Ky0Feb799sRiryIszF3i0Sj91CKe9+W\nui/X61baf+YeNE1/UfP82651QHgGkOer4DvKU7aWS0dxvfLLjjvSYaisH6V+zaPtneU8to8kW867\nJcm1rC/yjCKwVvunlN/GksnTjjjkhWhlWJ6aUs6pJVrbwcifq28z9dGnv+7ovdr/LW8xdXtSdKsT\ng8XfP8jzVRB1FMoiz2vTau9slV0jXeu5RHtndJ7rtSxf7vt0uJVuT7unFIszS+FK8bi3tsdtmXp5\n+ns9fdrWDyb/5aqTXrm+d272dJcs8oyizv0PkX3keTtj3e6Obl6kbu8d5HnnTNMvqN5RqMXRVOjJ\nzLf206r/HG2eL9HeuU5r0qxFoFLcCzeLQKW8N+779/t2z1I/izzLPZ7v3q0PuPbXHYmzyNJ2SvLC\nt7fI+LS1ve7sPcqotYO+PKtAp+n/1Dz/9pc+ITgI8rx7SsrWR5fRvI9IL42POi+Rso1SXq1ziKid\nx2jd6LzbbaFled/7ti6IKCK10vFPSYlSutGTWLxsM3mWVKqXViTN6LjlnKI0bSbQnvcmSon721Zu\np91zcySJ53zePQ/tKgAAAGAh8rx7ar1sr8HRezB72zuPRpaFnnPrT+HG6cHtefpIaS172tXx+85S\ntjZ9mbWH+vbCaFCF1u0t0r7N0+PbGqV9O2fpQWuPE43PW0vV2qnU7ihUS93eJj4bRLvnPYM8756S\nsq11GJKbP8pR2Z2C+do+1bHfUc5t8zwlYlxStqfgmzxqH8x74MYC9euzW1p8Z56aQLM0b4a/97NM\nrTSjDkCRPKPOQVmbZ9a+WWvzfLk2zuyg2/+Tdd5O4R4hbXvHTNPfJu0GgY+I1o1+cP2XQ7bsp9F8\na1+2nn31MrptT93ec9jWaUVHWc9RX2aFms1nQsraR71EW+cZyTE6RrTsz7fWEainfTN7j16Wnv8P\nK86TpukvXu904KoQed412UOve8iEldW9VKo2W9fba7eHvshy/bIdubYs+n76Mr+9d1TuWPvlVmRV\ni1Jr6c4oSvUCi46fYbfxw+llorTnHZ3DiBxrZX6blxtpKPsRlnVMI3V7ryDPu8b2srXivIRE7cvv\nK7v15Nx7PP3+Lsn5advt+t59bL9MIxGWeVvWSutm61v1vcBG8cfoiSKzc4rO3R/DH/vlo8uIIyb2\nn1m4N0jb3inT9LdoL09p+2G8xoczS8dm9Vr1s3To0ZRttG3rPI+kbe1yq72rr/wcmZX5SDzXEk5v\ndNybbs3q9Z7DyHaXZ+R/dPsZnab/49InA88Akefd8lNaB0Z4LnrSsCORau++7XY9jESWo9+y/vx9\n2rZH1tfBpiKfo9PMc3bMOSrD0RTteEq358dkVl5+3HK/5z2CPO+W99qnaY/+3O7tCBOlLiMZRiLt\naeM8R2o1/PUdlXR2/rX089HoOSeSZFbWeh09dm1/WVl2DaPURGrld922TSnPrmRZiaw5gATgPYI8\n7xZ/a8q5tKKlnja/nqgzk0xWbhmJKHu2Gz1G7YdCme9NUy/beOH4ebtcE1T04Gl7L6Z/GLV/6HRP\nZOefz5nte0SoNcHWzuMI1+8odBTkeY8gz7slG8e2R6bRr2Jfbl81WUTpy1rUmUWXmUSj88sYles5\nkXrrR0IceS6yWOftNJtvRY9WYNJ22T502g6l559e0nO7in3EWLbvSKI919CKUG9rfFqpP7KM6vvP\nAJ2G7hF+8twh0/Tzunzk2cORL4pa/VpkdolvS/8joLde7eW3UVC+n18FsRWnX9/zikTlX2Ws2egp\nJ/718eMyzQT19LTW+fhxnfcv/xzPmlSjCPVIWnkkHXw5Afu/e/Z/4evnTNOfP/+04Fkh8rxLSnun\nl2aPRMc/2Gu9WiTZ21moN207cm4ZvenakR8f/nxtpL3MPz09mXs9vXh9FLqPrmrpzZYspbYoHx7W\nh1p/uaqndWADm8K1Eax9GHaRqBVpmfcCLec0IlB/3SNCted+ObLPTSTQWnkE7Z73CH8xAACAQYg8\n75LoEWQRtSizlWry62rtnHbeT6N6hdE2z3PIjpvVyahF3Nso8+np6csrLlvPw0dbtkNPT8SZRZ4f\nP+6H1Yv9BK2HAAAgAElEQVQ6B5Vj9EaeHz5so89aWleK211rHY3se5GluGtctq20N+psHay2njjm\n3kCed4kfkq+H6INb+1Lw0pRbjlKwvt6RlG32BXMkJd3a/kj6Nk7bLqlaK8Nt6vZp8y3+9LnslHb+\nsYIYadOUthLrEefjY16nnE8kz+gVCbWck28fLdPo2pf3sJ62vX4notpnpuezI1eW7VOi09D9gTzv\nEt/T1lOLOKPlHlqRZCbLS0WdI1Fh3z6enk7JLRrxsUr9pycZMa7iLBHmsi76Et1GnrVes773bCbQ\njx/XJ5vY9sh37xaZRePNrueyTN+/X/fXK8+PH6VvvukTaKmfdTCKOhSV86u1hUbzdnqckc9PNj9S\nhjjvEeR5Z0zT36i9OM/58O2/4OP1o2lbX+bPsxV1XuPe1d5OST0/SvwPAnut63u6TdPqy7qnp0dJ\nD2H06dO2WYrWCq0IU1qlaZ9s4qVp5VxuXcnqluPbSLJI8ptv2hKV9ulcfz0+ArXneW7k2S/THtFF\nsmwJtLZt4aRp+nOSpHn++co+4FZAnnfHVzovQqulnaLUUyttG0WePbL159MzSMIo56Rqi/ysTPai\nX9b73rblhnwfhUq2DXSZP4XSjOTp2x5tSvbDh/2zM7/5Zo02s3ZMW788DDuTp488iyRbAi31fdto\nLQIt193bDmrPd5x1o3X73ugzE2j0WWrBE1buCeR5d5TOQhG9KaRMoBlZ9JnJ8JwOQ778XC4ZfdbO\nfRt1FnHuI89HPT09fC5b5qNbOLxEalGnFagUP+VE2gratpH6h1T7+jZtW7ax8rQCjURa6ke3t2QR\naDl+djuLP78Rlvo9UWZUVos8eyRZjz7hfkCed4ftKBQJydLzoe5Zl4miJyqN2kb9+Y1EnaMR6vnR\nZ1zXizQX6LK8mHCJNK1At9FnFIHZEYHsvZrRQ6mlvH3TR5FWnmU/rcizR54++pT2baGttk+pLc7d\nX6tboq2osqfOkcgz+zxmnwO4ZZDn3WF72mbUvkWyqDP6cHvRtaSYCVRJ/eh8R6/tSC/cXgHn4i1f\n1CVlu0Yzy3KRpo881+VHI9J3aecZaduZJxrwwL6WcwrOOpBnSan6nra98uyNPqW8J24m0XLdI4Mq\n2OuMI9Ja1NdbdsnI0y8jz3sCed4d53YO8tOeD3skxWy+p649B1XKldTJzjOjJcv+yHNt41zrxO2e\ny7yV6LL9VqqPj496eFg6DkW9aaWtHDNp2nRrLer86qtlH+/f79tIy8Otyz58W2wkT9/2GUWhUt4T\nN5OoPe8s6vRyrEWe2w5bu7WN5VqdVuQZldXOBe4F5PlqqIkx+zDLrY9STZfuMCRT35/LOR2hRuS6\n1l3EVxf2ttNQdC72fdi+j0WS23s+S9R5+tzm+ajT6Z0eHrbiLCIr82VYvSK52sOmy3mXaWm3tFGn\nlWftdpZI6GUfRYhZ9PnNN0v9qBduLXUrxbewZO2eGdv2zZeIPH0Z4nwtIM+7oxal9XwgW18g2Ta9\nUWYmVaktzOh8RiLto9FnbyeibH8lqtxGnz51u9Qp4nzS4+OTTqdHI9HTTqDSGnHaKNNLNCJr5yxR\nZyZPKY48I3l6gdro06dtR1K3Ui7NWmp2W7YX1dPOtOdGoudEnllduAcYEwoAAGAQIs87YZp+rrK2\nlSryZbWUbVR+cmW1CHM0bXtuBGgZSfsej3htqndN6e5TtyVlWKLLpb7vMHTS4+NJp9OSun18XPbn\ne9tmqdraQ6x9xFjaO+3tLdK+s1CUtrVp33JOPor0nYeiDkNlGo2D29PbNrv3007L+7z+LezUz69l\n9Q5GUdm5adv4nKbpz2qefyE4NtwSyPNuKDdQtzq/tNpXam0udl3W09auv3ba9kjnqEy256Vt9+2e\n230v7ZplX9sfGcs6f6vKMl3EuZY9Pj7odFqFKcXCHBFnlrIt8ow6Hdn9+H1J24EOorbP6FaV2u0q\nUYehWpunPb/sb7d2ztqWR3XrZUc+Y16StR+vHnrd3gPI825oZdizD6ev46dZm8stRJ4jUWe2j2xf\nttPQsm4rjmNtrX4YvjK/HSTBCnTtdXs6nXQ6PXwRo488rTBb4ozaOu1YuOW+TmnfWcgOkiBtbxUp\nbZ6241EWgZZ5KR+6r6fDUKu9c9t5KPoctGSVlbd+YEZlrf/Z+P9jC/K8B5Dn3RCJpfeD3PpCieoQ\neVaP8nn1OlD8ElXuH4K9zK+C9vd4nr5En8v0k06nd1+iz3IML9Ba79popKIiqXfv1rStHVQh6mlb\n9ulvVynH8GPreoFaedqyWuRp5elTtX70oX2HIRtpepH5tOw5Eu2NNu38yOcVed4DyPNuaEVURJ7x\nPrJ9jbR5Lutbqdu1nk3dajO/phL3Al166J42Eai0lWhGdPuGl6YdAH5UnlZg5Zx8ROtFGo1tW+tx\nW4s8s1tV7N9n/z6v02293bvX0d7ZEmhv9DnyeYVbBnneHUc/nPcYefr1LWrvSzvCXeVYE21y5M9t\nntv7RmvyfPrc5ln2f/qSql3bTrfHjVK2kVCshN69W6dWnn44P3+bSrT/aNSjlkCjgeFr7Z2+zdOL\nc3/NfvSm7fu9/Tv4eU/2IzJbrn2mWpFnTbrI9B5AnneL/1D2fPiyD2wkZFv+UpFn7VoijkSdrfX1\n4++j0dJJyG77pDK2rbSmepc6NvpclpfRipb07VL/tEnjluP6aZS2tdL0A8BHY+G20rY+Kozk5yUq\n7Z+mMjJIQn5vZ/x/H0Wdffd3ZutGBRp9nhSUI8t7Bnm+Cmq/eKN5XxZtm0WekUAVlHmpSrGUekXV\niv5Gth2R5XbbtVdtcAaf2zz3Eez272Ojz6UDkbQODl9SuFZk+ce0FnlaYUaDyNuh+GryLPu0aVsr\nz1YUKu3LWwKNri2KOH3kuX2/e4U5Uh7/Tddp7UeoX45ecC8gz7sj+8BFgoy28+tr+7KiVDLvJVkT\npz1uFHna8uj8a9Tk2i/LRYC1Y+6Pk6V74/FUF2GW7cq0CHO7/zL/SfZZj5lQsls8onTtiDx7j5G9\npFiWWdQZHcNKs/zwiGRZIn4bnbbaM9sdiWo/MrNp9lny+6uVwS2DPO+GTHjROg2ubx3XStTP97SD\n2uO20rU955SJt1U/2qYWlcbb1ToOrT1v17bPp6dHrbcalfk1An18fNLDw+NGTjHvvhzDHq/W5mk7\nC7UGku+RZ3SMWgQa9c714szStrVos/wfrkLd/o36ZZSty4TplzNxZhKufT6R5j2BPO+O2oc4+wD7\nX789v4ytJOWW/bykrsjTn/9o9Fnw++ypG+07+rI6uU4/fdvV2j73UemTSrrWbv/4qM9j29ai73da\nRytaxOijwtJJqOwvE2cUeUbXVab2GLYNNBKolaW0XY7q+gizbLOWPYUv+/fYRp1Pm/O39aLry2VX\nK8vE6ev3fBYjmcItgzzvht4PX/ar1+/D7yurH0n0FGzn9z3S7lkrP5fLRp1lXfnStSneRZSrfG0U\nmr9Xi0hXcS7TOAotkda7L+fw9LTWLzK0vWyzR5dF8txez+cjOpltU835y8qxbGtFmaVr95FnJMz1\n/cvf2+jv1irP1kVCbE2zz2tUhjzvEQaGBwAAGITI825o/ZLN6kXzviyr71O2Nurs6W0btXdm0V8r\nKhwh++XeTtnWz6UWEWzrbNs+t+Vru2cJLZd5G+FJebuqvqRvHzZtgvbWk1rU6QdE6I08fdq2FX3a\nKLIVafrU7fI+RBFnWVbQy9bPr+97Tm90OhJ9Rtmc6HMbbV87J7glkOfdUJOeL6+liWp1s7Ks3bPW\n5nmkvfNSqdvIOjU579fV2z23X9B+UIG47XDZZh031r7Xi0xLxyFpld8yv5WDlcfT07uwM09Lmn4Q\n+N42z6gN1N7G0tNzNpLlOu8HPfCv7XsXCzL639mWtXrg5uU18WVyrH3uou2Q5z2APO+GEqVkgst+\n2VpqEs3wEo0EqqDMS1XKhTUS2bXOtWf7XpH2RcXZvZ9xz9vPezg96umpRJ9rFPr4+LARqLRGktKT\n1gdsW7k8fN5uW/+IOP0PATufydNPvRTLteTCLNNHlYeEl/d1+0Nh/X+Nht2znX/2ctzWzWmLN5am\nr1cTY+tzizzvAeR5N/R++Hp+5Wbro+Vab1u7v1OjrDdlW+tp2qJ32x6R2p6yLbnHKVs7byPPtkBP\nmw49j49bua0p2jUiWzooPXzZ3grUyzPqXVtL25ZpJM+WRG00HNVZpmWoQj/ogRWoFer69yllMXH5\ntnrvtj0CbX22ap9hu69dTzG4QZDn3WAjzx5p9vzaHamfRZ2ZMK8RfR7leNS5jSrzc1wkeTLLT06g\nW7GuAi3lVqDL1LZ9PjwsEdnpZKPPJz08LOI5nR6+yLPML09qyeX55aoScdprK9OWPCM5Spk4ywD5\nT0acqwxt5Flk2TPge58cd1c5WOaPH312/Pqe18g5w0uCPO8GK89Hjf+ijT6c0Qfd0itQBWWZODN5\nRedwtNNQ7cunL+rsK9uWb9s7T5tlH3muAi2ytPO+Q9EajS6v5W+wRJtPOp1OnyVajrE8H/Th4UHL\n01lWidbk6efLNWXT7OVFauvb6NJ3CNqmZ/VleS/O9X33oszTtT3r1/3my7as5/OUSdQvl8+1ROR5\nHyDPu6En8vQfSk/vNlaQfrtW2rZXnD3tiUd/gfd0GOo5h/ULeo0qR+TvhbovLwJd7g8twiwb2A3t\ncR4+y7BEoqcvqVtJengokejjF3mWqXRK5VnrMFTmo+WaPG2kuMhyFWQmTp+OzUYQikWY/8+MRaSZ\nAKN12WfIlrU+d4+mDHneA8jzbvg8wnYYddY+mKqsy34x+/pZu2eUtm2Jsz+Si+kVY227XmFn55qL\nuDZw/LI+ag/9vIcvbaCSlWiR9yK9Mn36LEQ7LbI5fZHqOl32t5XoeqKtc7bTMp9HoU+hCLeR5lqv\niHDfm3YvUC/W6Dw/L+UXtN+yozwTaO3zZOd7pSohz/sAed4NXp6jEu0VZVTmpekFqqCsJVFVykp5\nxEjdVv3edSWyi463vx4vUB99tgV60nb8km00usp0L9BlH6fP84s01+jzpBLlbuV5+rKdj3a9MNf5\nrRh95LiXXlwna9Pcbruut8vRecX0pGu3+8zLsh+LvZ+vns9r+azDLYM874R5/hVJ0jT9rPrFWfsA\n+/KMmkCl7ba+bDRte4mUrWe0jdOu8wJdO/jk9UrdvUClvTTt/PYL3qdvy0AK5akvJ3OcItAiwWXg\n+SLMJercR5xbedr9rtNcnmWaSVRapLUKMO4QFIt2f7x9Wb8Ut9uM1+kR6CXkuUSc8/x3dZwjvDTI\n8+6w4uyRqBrLrS+AqN1TZp2XY0/0KbOclUW0OhCNbN8j7f15lZRovN+2bL00pX0nor1IpVWmRWon\nE00ut6msUlze/3X9w6ZsqWOjzFic9jrs+a/zVoxSkdm+w89+/T6ajB4ntj+Ofy97ZFgXbLYy+lHo\nlzOJRvX75Qn3AWPbAgAADELkeXd80vILtRZ1Kpk/GpnaSNLPy2wTdSC6ZJvnkcg022ak96wvq0Wq\nPtLcR6r7NtDSY/b0ZWqjzDVtW1ij0DUCXQdJKBFnWSc97VK2UZtnFnna8/bXsUaT67WUdX6c2f1A\nB3Fkuj9eFI0qYf9e9xFV7I0+s/mono80y9T2toV7AHneHS1xtuR4Tv1WhyEvySPi7E3h9tJqRx3v\nabuVX9SW2ydQqdZxyM/btO3+ZY+xSjOfLvvtl2dNXrE49eW9WuvtxRoJstWu2StEu+/OLTrKeuSY\nlbdepG3vCeR5d9jI04rUt4U+KP8AR2W9slUwL8XCHBHnJdo9R7bNRNpznk9Ofn0CXervJVrrRFTm\n29gfMtsodjst57e9VaUVda7n4SPROFrM1+0luX8Sit0m/tFxHrWdZLL0y5kks/maMMv0U/cVwMuD\nPO+OIs/oAzgqQ78uwm8TRZ8KptG6Uh6VZTKLzifj+dO2a6/atkCX+tuUrS+LIk9pTevWo8+y7ZMr\nj8W57r8vbWvPL1reS1Gb9yAqj6PN+vY9XC5dW1vOPkM9nz3/Y7e84F5AnnfHB+2jTv8hfFDfBzj6\n0Pv6Pi3rBapgWzutpW5HBHaUy6dt/fZjAt3fxhKVbXvqnlKp7iW6nntpS10FKq1CXevt07fb+RI1\ne/JeuPuycj6925qSpDznMrevRP/bfrklS78cidO+iDzvCeR5d3xUnq6tSbH3ldU/JfVObj6aSr1i\nWfd1DS6XtvV1zheoNoK0y/syG4na67ESLK9VsEWYXpbbdGwkUl83Pu9o+ei6nvJ23dH/o0yS2fxR\neVqJ2teHwfOFlwR53hnzPGuavqf9By+TqJLyI7K1645Gn76uKnX8+lFa+8sE6re9nkClOI1bBLkt\nW1OydnlN536prUyge8n688ui0KfNeWyvIzr/vGykA9BlI86enbWEaZfPledWpPP893ecH9wKyPMu\nyW5XiToNSe0PrpJ5/8qiT2n7xWyntk6pF5VFMrWMRhF+H72izM6pJ+XcEmh8rCVFu38PfCQqPWnt\neBQvL9v5bICPQJ/c8WoSXee3ctpez6hAyz766tW5TJo2W5ct93xuojL/WS2fY1K29wbyvEs+qt7r\nttXbtlekrfWRIDOp+nkNlp1LTaS97Z59UWku0OictvKJboHZr9PnqPVk0rnl77Fex9opaPujpjyF\nJU65S1vZbq8/vv2jP6W7Ld+L+vri3O2hUXYJeWZNLP4zC/cE8rxLbLtnLWVbPpAPyfoj8oyiT7lt\no6mflyLxxPU8WQq395t0JOKMytrLo71w99vtiSUqrVL8UqpVdOvfrKSBV3HWfuBEUai/RjXK82vo\nLa9xmdtWpFyQ0fIRefrPZNTTlsHg7w3keZd80D51m0lUQVm54f6IPKMvjbJcS9mOCsjTSsH6Oq26\ntWP2tnvaOkcF2pOuVbje14l6qu7bTf3fScm8L5Piv4GX5bZuS3BHBRinkHe1RvZYWb6kPL00kee9\nwti2AAAAgxB53iFLj9ufUdxxyP6atQ9WrqVubSSqSr3yOrl5u43MNOs05JdbEWC0fS+9nYai9T3R\ncTuKrD1oukYWYeZ19hFi3ObpswNZJPrk9rVd3tb15xTVjc57vH17f99p7/a1eteKPLMRwOz8J83z\nP9B5DXArIM+7JWr3rKVt/diomRhVWZfVk1mOpq32s0xqR4WZ0WpTPdLe6cv289sBCfJ68bzMPuoi\n9unTU1I5T5X2yHEv6OX8spR5JNYjZG2soxwR6BF5PiZlvpmFzkL3CvK8W7Iet9GHMos6o7ZPJXXP\nkWdNlEc7DZ3DpTsM1faxrbsdJH70B8J6b6jUjmizCHC9rSVq81RQ9qRWR6D+SHNzJuq7/u11n0+2\nr3Pk6X+o9nyWaO+8Z5Dn3fKNFnlGHYf8B7qVtr2kPM/pNKRg2TOaA+2VXHSM2rm2OwzlAo3EV0//\nRvSkdU3tL/vNhtvb1uvdX/vcjmN7DF+aHoGOyjP7jGTNKsjznkGed012y0oUeXqB+qizLEv1L4SW\nPGvCbImyJo2jEWlPqrb3fMbbPLNjrqML9WzTruNldbSd1e1V9eg0P36b9t/6stFmRu1/8lLytJ9J\n//n8JOnrcy8CXgDkeacsnYZ+SvvoM4o8vSBtu6eXqFT/QhiVZ/TFEy2PRIFH6NnPddo84+Pb9lAF\nt7SMpIb3RN65hFDHfFb7m8Y7WqLioz+QvuxlYNtry7M1MMInzfPvGDhfuBWQ513T0+4p9adtS/1a\nh6JMnpGAFaxXsjz2JXuca3Uayur3RZa1wRHa59fHswRyB1k7VfWeZC1joM51UZ0RedYEajsMZZmh\nTyJle78gz7vmo1aBZtGnTNlJcep29FaV7AukZ6pguUdYl6aVNuwRqExZf8o222Z96sklfkhEbc+3\nwP58xsR5JBPRqntpedoOQ1HkSXvnawB53jFL6vY7qkefUixMm7qNUk2t6FPBfM9UwfJRcbZCtXbE\nl9c7mj6MxDiW+t3eA6lEpnYbO1Vy3s/9Y2Qkrdzb7quOeukRO9eNytMv+x+ukTBtyvYfGroKuB2Q\n591jI89Pkt4pTttGkWcUfcqV98rTpoZVmSpYrq2L6P0i7U3tZfsbbe88QnuftUHjL3O8nrqjqet4\nf/tOQJfpMNVHto9z5Rm9avd5krJ9DSDPu+eD6m2f0vrBLb1uT+5V+xUdCdS3ocptp8q0le48r6PM\neds9f4ehvnMsnYvW9ZfpTRtxCVFFUfRLnk/hUvJsCdSnbeOok4df3zfI81WQRZ8+Ksw6DtnoU9oK\nM7qlxe+j0CPPVpR5DXHa7Y8c4/odhkYj2rF7PHuPeZQ4Um4f71o/lDIuKc9MoLZTkJLy8jn9yZGL\ngBuBgeEBAAAGIfK8c5ZOQ+8V97j1bZ5Z2tbfZhJFm4Uoci1cM/K8XqR0mWOd39t27NzWafR8zu3g\n7+de76U7/7Tqnds5KOOcyLMVffpRhGpp26W9c57/sUNXAbcB8nwV+FtWbOpW2n5wo05DNj0rt5y1\n63ix+jbQS8vz0l+kR2/9sMuq7OMSqduR9VtqqdPR+0njZ2f2XtuoCK/ZK3hEnrX/29bnwX4u7K1j\n5cctnYVeA8jz1fBB0rdUjzxb0adv8+y9ZUXafnkc7XUbLV8rCin7Hol0RzoL9R73+SOwfBSfa9/O\n0jrOtY9/CXn6aDQTp49KrTg/Sfrx0YuAGwF5vgKW1O2D1gj0nZYPaGnSLsLsuWVF2kqzJtBHcwwf\ndWZfOjZSteXqXL40owId2V+tt+3Rc4rq2WnvMa7BJW5nuSbRuUTitPPZj0Qvx1bk6VO2/+Thq4Db\nAHm+GspoJb7tU9p/gLPo00aeXqARPsK0UWcWgcotv7Q8pfPFdqk20tH2zp79H4lse2696Tnf2nn1\nll+SXnleImUbffbKZ5NbVF4DyPOVsESfX0n6SsuftbR7SnGnhSJML1JpH43WvjCydJbMPuxtMgrq\nKyjPlmtEjXnnfiE/x5f6c9yy8dxR4K1FndJ58mwJ1Iszk+dHzfPvPOsq4DZAnq8KO2BCeUmLtLLU\nbRRhnlxZT7tnNHBC6xe83c6XZ8sZl4hoLtFT1K8/Gq1dSz7XivzO6RT1ElFnbTn7IWjLMoFmaVui\nztcG8nx1lHbPEn1K+yH7aqlbmXJpm84ty/bLouDHz611Guqdj5YjLhm5HU3LPlfq1teL0rOXOt9L\np2xfkl55nht1Rmlb+2P2m3MuAm4I5PmKWFK3Jy2pWxt5ls5DUfQZydO2efbc8ym3vqfd03c2UjAf\nLUe8ZCR3KS7Z3tm7z0ttc84+bqmt0877qLMVeWZpW9vL9oPm+XedcyFwQyDPV8eTtr1upbX9s0i0\n1utW2ovT38ISpW29LKPynnZPX8f3zs2uuYdzvqijyKt27EulbmvH6OXagrr1HyaXiDpr8ozEGUWe\nRJ2vCeT5ylhvWykdh6T1wxtFnj4ClfLI06dvfS/co6nbVrtnVlbw25/LURncukRess3xpd6bXnHa\n+SzStPWy21PsPZ2Sbesk6nxdMLYtAADAIESer5JHLb36yp/XpmwftL+3M4o8a+2eWWR5rXbPrGxk\n/a1yqftGL81oB6Rb5VIp2yxtm6VsbX+Dj6KX7esDeb5aPmhJ3Urb+z5r93x6eUbtnmV91uYpxeK0\n7Z6ljsz2VqRRGraVmr106vYeONKJ6BLH1DMe7xx6UrR+uSbKqKxXnr96xnXALYI8XyFLu+ek9dfu\ne+2jz155yk1926cla/P00acUS/Pc6LMl0CKaS1ITl4/eVKnb2ldPvV6xvWTb53NxtK2zTLM2Ty9N\nL05/j/USdc7zP3P0QuBGQZ6vmiLPr7TI8522qdsHLR/ymjx7et3KLWfizKLLS0WfpU65Dt8J6dLi\nfG5u7RaQW6X1w6slzjJtRZ5ZRyEvTwaBf40gz1fKGn1Ki0S9PK04bfQpxQMoWHFmkWehJs5zo8+R\n9GwUVd8Kb1lulkv/qBkRp13uaeMsy9mgCFac5ZFjHzTPv+fIhcCNgzzfBN9oTd16gZYPvBVMdBuL\nlEef0v5Xei36lDn2aPRZK4c9tyLpa0b+tetrRZx2vvb/HKVtrTw/uVeRJ1HnawV5vmLmeZYk0/4Z\nydOnZqWtTB91+ehTZj4aQMF+sWWDJFxToLeW4r0VAV6Lc97vS4qzTLOXXR91vCsS/agyIMI8/96R\ni4E7Anm+GXqjT2kv01b0KcWRp9SXtvXT3vZOItDXg+9Y1Vu/Z33tf86Xtf6fM3naiPODpB91Xgfc\nK8jzDbC2f36jdvSpz/OfzB6ynre16LPW87asr0WfLyXQ0SjoOaLUW4uEj9JzHedG2Fmbpp0/J/Js\ntXMuUec8/3NnXgfcOsjzjbAKtESfkTyz3rat+z5ltnsK1kX3d2bifEmB1r7cjwjs0vs7wrU65NT2\n+VLX3SPOaL5HnJk8fTvnB83z7z73QuAOQJ5vjtL26eVpBedvX8nu+7TzUaq2FXn2pG2jNs9MoKeg\nvIfRlOERnuMY55KJrSW8l4yMs0h1pI2zTLO2TS/PWtT59ZGLgDsEeb4htoMnlPTtR63SsXLMBCo3\nbXUckuIOQ3a5J+psLdtye349gyb0pBJ9nUiG50rkkpI6N/q75rVc4ofEOR2FnoLyWltnJs8o6vx9\noxcCdwq9LQAAAAYh8nyT1DoOSWs0FkWfCqYWe/9nre3TLl86+pSrU2vDbXGLnXVq6VUl6659Dll0\n3oqk1aiTbdOzrjddW6bRy7Zxysz7qPOjpJ8MXAfcO8jzjbFP3dbkaXvdRu2c/kvPp3F9j9tSJlPW\nI1CpLtSjbZ0tjrYB9tYZqVfb9trp4nPo7WB0iWNEZSNp2x55ZinbrzXPv//M64B7Anm+Wbw8W71t\nP7l10l6oPbeu9EaetSjz0h2GIkbFOdLG+FIdcEaln81n+zu3nXWUmjSz+RFxZvL04vxGPDXl7YE8\n38m2zVIAABEQSURBVCD7+z5tj1uZ+SO9bltEEeVRgUayHOkwFNHbiSjaJls+yqXSsM+Zen6JjkI9\nkedoutaLM5JnGRDha83zv3TgWuCeQZ5vlK1AbepW2kec2ZddTZ627bPV/tkSaKkjxUKtSbLV5hlF\n3BFHpGq3a627RjR6iWh4NPpsHbesV6NOVL9nXW87ZzRve9ZaaT5pzbwUedoBEf5A53XAawJ5vnk+\nqC5Pv5x1uGmlcLO07Ujb58mUe2G+RNtntO5cCfauz6aj5ztaJxJfJlBfL9vXEXpTtrVIU67sMZi3\nbZz6PL+Kk3Tt2wV5vmHi9K20SqgWndW+FGvrfNq2lB2JPGsSlc4XaS3a7JFIb7R4JHobpUdw12jv\nvFTq2e6rVtZq5xyRp7+nU9qma3+ief5XjlwIvAKQ5xtnn76VcnmWNK4lij57viiP3LaSRZ6t9s9T\nsj7CCrOWxhyVSU1U2Tat8hpRRHqJ/dXKeiLWQu+PipF1I+2bcmVZ1GkfcC1t2zn/xY5rgNcK8oTP\nrI9R6k/TRu2cI3iBnlx5T+R5qfZPKb4my1G51SRzRKK1qSrnMtqu2SPH1o+EjCMp25HIsyZOL89W\nura8pFWe84Hzh9cE8gQXfUr7ez8jydQi0Gg5I+o4ZOdbkWeWui3ncG47aG8kmpWdm7pV5/ZH9t2b\ndm3JsXbuhaPXkEk2Ss9G85k8o3IvTt+zVlo+Iz/WPP+h8UuBVwXyBElWoNIyuHUWcSpY7l2XkQnU\nR7q1EYhqUWcUcWZC9bJsCSWSSG9Ueemo0+/3aMr2GhHnkUiztY+jUWckz5Y4P3yu/7Xm+Q+eeyHw\nCmBsWwAAgEGIPOEL87y042w7EPW2EUYRx0gUWrvvMyrL0ra9gyT4aNRfT0TWWciuqy1L++gwqtuz\nvpdWFNp7Tr1R7KVSzX5/WVk0Pxp5RlGnvZ/zg0qTxjz/y4evBF4XyBMSPmn5wrimMGvblvmaFFtt\nn3a7o+2fvW2eLXFmy73p1mh9rcNQTXaj55NdRylTpdxzNK1bS9lmAu2RZ02cn7SK84eN84O3BvKE\nHdvB4yN5Kpi/BFnbZzlWT/TZijyzHrdKykfaPP261nJP9FmOcST6HGn7rJ2XgmVbViuPzqmX3qjz\niDwfzbpMnKUH+o80z3944LzhLYA8IWTbA7clk0uSyXlEnNHtK0cjz8JoBBotK1ifzR+NPv3+su1G\nzqm2T1XKLecINYs6e+Tplx9deU2cZSCEf61xfvAWQZ6QsgrU9r4tXEumtci2N21bu1XlSOQZXedT\nsI1dd6no05YfSdv2iPk5os5rRJzRfCtVawUqxeIs6dqfcEsKpCBPqLKPQKXzpNn6Eo2iTb+uR6Jl\ne+l2ok9b1jN/jbRtr0AVbOP3rUp5RCbcjCzitPOj8rRtnFIuzq8RJ1RBntCkHoGeQ/Tl6MUWRaKj\nPW97ok8FZa3os5yzX59JVEnZNdO2LYEWvKyj/cjVr5VH1ISZrWtFnJE8a6naSJ5enN8gTmiCPGGA\nRy0CLVxSopE4/TGiARR607YjAyXUzjMSbrTel2V1vNxqUeGRtG20fe28onl/nT3lR+lN2XpByi1n\nUWeRp3/EGD1rYQzkCV2s0ae/hWWELIqo1c8i0WjaSut6aT4G21068mylQntE1yvO3mizdU7ZOktW\nbte16Ik2/XJNnlHKtvakFGkrzmXM2nn+I53nD28Z5AndrAL92Kzbj48cbLkXpy1riVLaizNrA/Xy\ntMf1IjwaeUbStMt2/pKRZy1dG4myts7iU77Zuh5qsvTLvdKMBGrbN4s8y7M5P0j6Vc3zvz147vBW\nQZ4wxPkCrUWf9guvlcLtiUDt8kjqNhp9KLuOSKheuD6K9HWi+UtGnrV0rj8nBdvYdZbaulFq6doo\nLevn7XIr4vRPSSmDvSNO6Ad5wjBjAm1FFVlZVMeKblSedj6KQv36npc/v5pEo3qteb/N0cizVqbG\n8pEIs/WDo3ddFFmWaW3eS7PMl4izRJvSEnH+SPP871TODWAPA8MDAAAMQuQJh8ijzyiVZsvfBfX9\ntv4l7dshW22fUnwbS5bClVvujTxbKVtbpqSs1e55Tto22l5u3i9H6VhfX26dpZVJaGUfsmxFLcr0\nZa1HjJVn1/5Y8/xHG+cLsAd5wmG2zwAtEs3apyKsTFt17TYPwbwXjbQ+FaaU10RZynqTMVmbZybR\nqH5Pu6eS8p60rTrWZ/XtOin+29TW1RiRZ0/Ktiz3tHOWzkE/kiTN8783eO4AC8gTzmL7GLOPir/k\nFJT5KLQ3qojEWebL1EZCI/LMOhXVIs9y/Gh9FIX6+q15ufmWSEt9BesVzEf1lSzbbRSsGyGLLrP5\n2v+Hl2etnbO0cf6xM84dAHnChdjeB/oTxRFoFl0WkdbE6/fVij5l5q1UI3nWbmFpybMWbfbIrzUv\nV+731bruVkTq62fLWdkRLi3OSJ6ZOH9V8/wnLnAN8NZBnnAxVoGWkYgiYUZSfefW1b5As7bPSKI+\nEs3kqaBsVJ4tiT4F9XuiTgXltamC+ajspcSpYD/+7177+0ev6PFiXpzLyEHz/B9c6BrgrYM84aJs\n07glAs0iT1/u2z+jaNS+ailcKY7Mri3PSIaZWP32UYSqoLw29WUtiUrbv8lLizMqa4nTRp6ROL+m\nbRMuDvKEq7AdTD6SZ+0L0UuzVf8hmJe27aA2CrXzJ1f/UvKM1ikoz2R6NPq071lNmH7e16uVHSUT\np51v/cCKxFkb5P1rzfMfv9D5A6wgT7ga28eZ9cjTlr8zU19fQVlPR6KaQBWUnSPPKPpUUO73kUWZ\nR6LPlkQjUXqisiO0Is4sAm1J045TW+RZ2je/1jz/yQudP8AW5AlXZRXoB22/6OyXoH29N/M2An1n\nysv2D8H2WeQp7aVp56P158rTr1dSXos+/TY9U1/WE3WOlI0QRa01cbbkmd2/KW3F+WPN83945rkD\n5CBPuDr7nrjS/oswkqaPPG35w+flTKAy87Y8E6i0lVaPQDW4Xo3yaD7appa6jbbxZX5+tKyHmjTt\nfCTOKDvRug1FWsX5Q83zf3LwvAH6QJ7wLGw7EkmLRKPIM0vn2jSuF2UkUGmfsi2SiVK40r4d9KS6\nRFVZd015+jK5+ln5EXnWyj0tYdrlmiizsmjgAztO7df0poVng7FtAQAABiHyhGdlG4F+o3004VO4\nURtnSdm2olWbzow6Dvk2zyjarEWgCspq65SsU2XebxPV9fWibfx8tJyVHSGLOO18LcosZVEbp+8Y\ntDQFzPOfutC5A7RBnvAibAeWL7cWfEvx+KT6PPWp25pEpa08ayItwrBp3EyiJ7NOyiWYCbW1jS4w\nPTKvjvIWNWHa5ZowH5XL85OZftDavvmfHTxfgOMgT3gxtlHo11q+GItASwRqe+f6ts6aRMs2J7eN\nF6kVoZdmLeosx6jJ8qg4LyHPrHw02qyti9o4ffloxOkHPbDi9Pdv/seVcwO4LsgTXpzt7Sw2Ai0C\nlbYytVGol2h5yczbL2YbbdqptIoyk6gXqbSVXyt1qwPrZJajOnL1svJrp2xbwrTztXS7l6eNNssw\nez/SPP/nZ5wrwPkgT7gJ9sP6lSjjq881ijitQItEy7TIMhpYwYqzJk8rzkiatkzaCrAn+jwqzmtE\nnb2ytO9RRo8wy7QmzmjQAxtt/unOcwa4LsgTboptFFraQ6VFoj6dG0WhfpCEB/VFoDLzUbTp2zxP\n5hi90acq666Ztq3NZ/Vb9HYIsmWtaLM26MGsef6vBs4P4LogT7g59m2h0jYS/aTlXzeSaE2eXqSP\nrn4mz6izUBFNJNWaSJWUZxGqkvKaOK8pzUJPZyA/bUWbNtKUSk/aef4vDpwfwHVBnnCzFIlK9mHb\nH7W0iZYvWS/R8kBraTsKUUuiMtvaaSRP20P3SNo22zZaf2SqYDkS5BFpFlqp2TKtdQwqEWf0+DBp\nnv/MGecHcF2QJ9wF22i0CLREoUWi77S2gUpxe2j2kvIORVakl5Jnb+pWSfm5Eeclok1bVhNmmfdD\n7O07AxFlwr2APOGu2D7q7IO2ErUClbYybUWh0j7ifDDzj8G89HzyvGTUee00rX+V8tpABz/SPP+3\nB84L4GVAnnB37HvmftAi0K+0pnGlVZ5eopFIpa0wI5F6gdpt7kmerfJCrVNQWe5p05SZtyna8vST\n/65xHgC3B2PbAgAADELkCXfLNgIt9wKWCFRao1AfgdooMhpUwaZra1NpjUJ7IlA1ls+JPGvz0XIv\n5w54kD2o+sea5//h4DkBvDzIE+6evUS/+bzmW1oFmkm0TKW8PbQlT5vW9TK0t7UoWF8T59GpZaSz\nUK1DUJlv3bsZjUUrbe/X/F+S4wPcD8gTXg321hapyPSkWKLvtZdnrUOR70AUyTOSphWqVBdmK/rs\nmfr6rTJLT6egSJ7Rk09slCnN83/fODbAfYE84dWyH2zBdyqyUai0j0ZbApW20szma6JspXRHpkqW\ns7JCa7SgnltPrDTLUHr/c+WYAPcN8oRXz36whQdtRepvbWlJNJKnn0YSlfrEea10bYvWPZvR48K2\nUeY8/08HjgtwfyBPeFPso9H3qncwqklUitO5WfQpxfIcTdtG5X5eSXnUrmnLs9tPbJQpbQds/7Hm\n+X9L9gvwOkGe8CbZR6MFm9a1Is3uDa2ldXsjz5pAe6Z+vofsns1S7h9KXmT58fP6H2ue/+zgMQFe\nD8gT3jy+o5FUhPpedZFK/b1zozbS0cizV5q+rLcjkI80V2HO8/8aHAfg7YI8AQLiyNT22u3pZOQl\nKtU7FF1Knp7a7SWlzHb6+ajlaSb/e8e+Ad4myBOgQR6ZSnHv3dr4uaOdiuTK7HIPWfulfW7mR83z\nn+vcHwBIyBPgEJFQJZvu/Up5mjdqGx1p+/TzUj76j32CyQfN81/ovEIAqIE8AS5IJlXLNP2M1mg1\nS/P2RJrb9sl5/sVzTh0ABmBgeAAAgEGIPAGemXn+lZc+BQA4EyJPAACAQZAnAADAIMgTAABgEOQJ\nAAAwCPIEAAAYBHkCAAAMgjwBAAAGQZ4AAACDIE8AAIBBkCcAAMAgyBMAAGAQ5AkAADAI8gQAABgE\neQIAAAyCPAEAAAZBngAAAIMgTwAAgEGQJwAAwCDIEwAAYBDkCQAAMAjyBAAAGAR5AgAADII8AQAA\nBkGeAAAAgyBPAACAQZAnAADAIMgTAABgEOQJAAAwCPIEAAAYBHkCAAAMgjwBAAAGQZ4AAACDIE8A\nAIBBkCcAAMAgyBMAAGAQ5AkAADAI8gQAABgEeQIAAAyCPAEAAAZBngAAAIMgTwAAgEGQJwAAwCDI\nEwAAYBDkCQAAMAjyBAAAGAR5AgAADII8AQAABkGeAAAAgyBPAACAQZAnAADAIMgTAABgEOQJAAAw\nCPIEAAAYBHkCAAAMgjwBAAAGQZ4AAACDIE8AAIBBkCcAAMAgyBMAAGAQ5AkAADAI8gQAABgEeQIA\nAAyCPAEAAAZBngAAAIMgTwAAgEGQJwAAwCDIEwAAYBDkCQAAMAjyBAAAGAR5AgAADII8AQAABkGe\nAAAAgyBPAACAQZAnAADAIMgTAABgEOQJAAAwCPIEAAAYBHkCAAAMgjwBAAAGQZ4AAACDIE8AAIBB\nkCcAAMAgyBMAAGAQ5AkAADAI8gQAABgEeQIAAAyCPAEAAAZBngAAAIO8v8ROvv/9719iNwAAAHcB\nkScAAMAgyBMAAGAQ5AkAADAI8gQAABgEeQIAAAxyenp6eulzAAAAuCuIPAEAAAZBngAAAIMgTwAA\ngEGQJwAAwCDIEwAAYBDkCQAAMAjyBAAAGAR5AgAADII8AQAABkGeAAAAg/z/AoIWKJszDagAAAAA\nSUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0xcbb4a90>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "img = run(w, h)\n",
    "fig, ax = plt.subplots(1, 1, figsize=(10, 10))\n",
    "ax.imshow(img)\n",
    "ax.set_axis_off()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "8.33 ms ± 14.4 µs per loop (mean ± std. dev. of 7 runs,\n",
      "    100 loops each)\n"
     ]
    }
   ],
   "source": [
    "%timeit run(w, h)"
   ]
  }
 ],
 "metadata": {},
 "nbformat": 4,
 "nbformat_minor": 2
}
